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Abstract 

We have extended our method of grouping of Feynman diagrams (GFD theory) 
to study the transverse and longitudinal Greens functions Gj_(k) and Gii(k) in (p 4 
model below the critical point (T < T c ) in presence of an infinitesimal external held. 
Our method allows a qualitative analysis not cutting the perturbation series. We 
have shown that the critical behavior of the Greens functions is consistent with a 
general scaling hypothesis, where the same critical exponents, found within the GFD 
theory, are valid both at T < T c and T > T c . The long-wave limit k — > has been 
studied at T < T c , showing that Gj_(k) ~ ak~ x± and Gn(k) ~ &fc _A n with exponents 
d/2 < \± < 2 and An = 2X± — d is the physical solution of our equations at the spatial 
dimensionality 2 < d < 4, which coincides with the asymptotic solution at T — > T c 
as well as with a non-perturbative renormalization group (RG) analysis provided in 
our paper. It is confirmed also by Monte Carlo simulation. The exponents, as well 
as the ratio bM 2 /a 2 (where M is magnetization) are universal. The results of the 
perturbative RG method are reproduced by formally setting Aj_ = 2. Nevertheless, 
we disprove the conventional statement that A^ = 2 is the exact result. 

1 Introduction 

Phase transitions and critical phenomena is a widely investigated field in physics and 
natural sciences E3 II]- The current paper is devoted to further development of 
our original diagrammatic method introduced in [5], to study the </? 4 phase transition 
model below the critical point. Our approach is based on a suitable grouping of Feynman 
diagrams, therefore we shall call it the GFD theory. In distinction to the conventional 
perturbative renormalization group (RG) method 011]) it allows a qualitative analysis 
near and at the critical point, not cutting the perturbation series. In such a way, we have 
found the set of possible values for exact critical exponents [5] in two and three dimensions 
in agreement with the known exact results for the two-dimensional Ising model |7j . A 
good agreement with some Monte Carlo (MC) data [S] |5] and experiments |1U| has been 
found in 0, as well. The disagreement with the conventionally accepted RG values of the 
critical exponents has been widely discussed in jlll I12j . providing arguments that very 
sensitive numerical tests confirm our theoretical predictions. 
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The ip model exhibits a nontrivial behavior in close vicinity, as well as below the 
critical temperature T c , if the order parameter is an n-component vector with n > 1. The 
related long-wave divergence of the longitudinal and transverse correlation functions (in 
Fourier representation) at T < T c has been studied in |131ll4j based on the hydrodynamical 
(Gaussian) approximation. Essentially the same problem has been studied before in ^Hj 
in terms of the Gaussian spin-wave theory ^H] . Later perturbative renormalization group 
(RG) studies [IZlIlHlIinil^lEllESlEHlEllES] support the Gaussian approximation. The 
RG method is claimed to be asymptotically exact. However, we disprove this statement 
by finding out the errors as well as unjustified assertions and assumptions in the papers 
where these claims have been established fSec, l9,l|) . Our analysis predicts a non-Gaussian 
behavior, and we show by general physical arguments that it must be the true behavior 
to coincide with the know rigorous results for the classical XY model. This prediction is 
strongly supported by a Monte Carlo test fSec. ITU|) . 

2 Primary equations 

We consider a <p> model with the Hamiltonian 

H/T = J (rV(x) + c(V^(x)) 2 + V(x) - Mx)) dx , (1) 

where the order parameter (p(x) is an n-component vector with components </?j(x), de- 
pending on the coordinate x, T is the temperature, h is an external field. The same model, 
but without the external field h, has been discussed in [3], representing the ip 4 term as 

(/? 2 (xi)u(xi - x 2 )<£ 2 (x 2 )<ixicbc2 (2) 
= V~ l ^(ki)(/5i(k 2 )n kl+k2 v3 j (k 3 )v? i (-ki - k 2 - k 3 ) , 

ij,ki,k 2 ,k 3 

where in our special case of (^Q) we have u(x) = 5(x) and = u. This is obtained by 
using the Fourier representation y>j(x) = V~ l l 2 J2k<A Pifo) e lkx > where V = L d is the 
volume of the system and d is the spatial dimensionality. Like in here we suppose that 
the field ipi(x) does not contain the Fourier components (pi(k) with k > A. At h = 0, 
the model undergoes the second-order phase transition with a spontaneous long-range 
ordering. Besides, all the directions of ordering are equally probable. To remove this 
degeneracy, we consider the thermodynamic limit at an infinitesimal external field, i. e., 
lim lim where h =\ h |. In this case the magnetization M = {(p) is oriented just along 

h— >0 L—*oo 

the external field. We consider also a model with Hamiltonian 

H/T = J (?V(x) - 6 ■ <^ 2 (x) + c(V(^(x)) 2 + u^ 4 (x)) dx. . (3) 
In the limit lim lim this model is equivalent to the original one at lim lim in the sense 

<5— >0 L— >oo h^O L— ioo 

that the magnetization is parallel to certain axis labeled by % = 1. Some degeneracy still 
is present in Q, since two opposite ordering directions are equivalent, but this peculiarity 
does not play any role if we consider, e. g., the Greens function Gj(x) = (y?j(0)(/?i(x)) . In 
the Fourier representation, the correlation function Gj(k) is defined by 

(^(k)^-(-k))=^G i (k). (4) 
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Hamiltonian (J2J) is more suitable for our analysis than lfl|). since the equations derived 
in |5j can be easily generalized to include the symmetry breaking term 5 ■ ip\ (x). It is 
incorporated in the Gaussian part of the Hamiltonian 

H o/ T = E ( r o " <* • ki + * 2 ) | ¥><(k) | 2 . (5) 

As a result, the Dyson equation in [H] becomes 

' -r -d-S iA + ck 2 -^^\ + ^(k) (6) 



2Gi(k) u aGi(k) 

where D{G) is a quantity, the diagram expansion of which contains all skeleton di- 
agrams (i. e., those connected diagrams without outer lines containing no parts like 

— ^ — ), constructed of the fourth order vertices > <C . The solid coupling 

lines in the diagrams are related to the correlation function Gj(k), but the dashed lines to 
— y _1 Uk = — V^ 1 J u(x)e~ ikx <ix. Here the notation = uu^ is used for a generalization, 
while the actual case of interest is = 1. Any two solid lines connected to the same 
kink (node) have the same index i. According to the definition, Eq. © is exact. It is 
ensured including the remainder term i?j(k) which does not contribute to the perturbation 
expansion in u power series. Quantity D(G) is given by 

D(G)=D*(G,l)+ O-O (7) 

where D*(G,Q is the solution of the differential equation 

D*(G, C) = -~£ Ml - 2E(q, C)] " C^*(G, C) (8) 

with the boundary condition 

D*(G,0) = -^ln[l-2S(°)(q)], (9) 
z q 

where 

S(°)(q) = -2n q y- 1 £ G i (k)G i (q - k) . (10) 
Here S(q, C) is a quantity having the diagram expansion 

£(q, = - q -0 : " + C - q -/!V- + C 2 I q -/Tl>^ (ID 
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including all diagrams of this kind which cannot be split in two as follows t $&T xTx ^i • 

These are skeleton diagrams with respect to both the solid and the waved lines (i. e., do 
not contain parts like ^ and/or vyxp^xrv ) with factors 

kr^rx£ =-n k y- 1 /[l-2S(k,C)] (12) 
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corresponding to the waved lines, and factor —V~ 1 Uq corresponding to a pair of fixed 
(formally considered as nonequivalent) broken dashed lines in Quantity X(q, £) is 

defined by converging sum and integrals (cf. Sec. 4.7 and Appendix A in [SJ), i- e., 

u~P u-P 

S(q,C) = S (0) (q)+ J e-^dh J e~^B{^ Q, tfa) dt 2 , (13) 



B ^c,t) = Efen2 s(n) (q>0, (14) 

n=l \ 

where S^ n ^(q, £) represents the sum of diagrams of the n-th order in (jlljl . and p is a 
constant < p < 1/2. Note that the zeroth-order term is given by Eq. (|1U[). 

Based on these equations of the GFD theory, we have derived the possible values 
of the exact critical exponents 7 and v describing the divergence of susceptibility, i. e. 
X oc (T — T c )~ 7 , and correlation length, i. e. £ oc (T — T c )~ v , when approaching the 
critical point T = T C from higher temperatures. These values at the spatial dimensionality 
d = 2, 3 and the dimensionality of the order parameter n = 1,2,3,... (only the case n = 1 
is meaningful at d = 2) are jS] 

d + 2j + 4m 

7 d(l + m + j) -2j { ' 

■ - 2 ( 1 + ™)+^ (16) 



d(l + m + j)-2j 



where m may have values m = 1,2,3, .. ., and j may have values j = —m, —m + 1, 
— m + 2, ... A general hypothesis relating the values of m and j to different models, as well 
as corrections to scaling for different physical quantities and several numerical tests have 
been discussed in |111 112j . Here we only mention that m = 3 and j = holds at n = 1 to 
coincide with the known exact results for 2D Ising model. 

Since our equations contain the diagram expansion in terms of the true correlation 
function Grj(k) instead of the Gaussian one, they allow an analytic continuation from the 
region tq > 0, where they have an obvious physical solution [5], to arbitrary ro value. 
One has to start with a finite volume, considering the thermodynamic limit afterwards. 
In this paper we have extended our analysis to include the region of negative ro values 
below the critical point and to study the transverse and longitudinal fluctuations of the 
order parameter field in presence of an infinitesimal external field. 



3 Alternative formulation of the basic equations 

An alternative formulation of our diagrammatic equations can be helpful to prove some 
basic properties of the solution (see the end of Appendix). Namely, a resummation of the 
self-energy diagrams contained in -Rj(k) = —dD*(G, l)/9Gj(k) can be used instead of JHJ). 

As discussed in j3], these are skeleton diagrams like X with two outer solid lines, 

obtained by breaking one solid line with wave vector k in the coupled skeleton diagrams of 
D*(G, 1). According to our notation, the factors of the lines marked by crosses are omitted, 
and "coupled" means that the diagram does not contain outer lines. The summation over 
the linear chains of blocks ^ contained in the self-energy diagrams of i?i(k) 
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can be performed, as it has been done in [5] with the diagrams of S(k, Q. It yields the 
expansion of -Rj(k) represented by skeleton diagrams where the true correlation function 
Gj(k) is related to the solid lines and the dashed lines are replaced by the waved lines. 

Like in these diagrams do not contain parts ^ and Kf\r^f\nj . According 

to (JZJ), the diagram O — O i s n °t included in D*(G, 1), and it also cannot be involved 
in the actual grouping of diagrams since the extension of dashed line by adding the blocks 

^ would yield non-skeleton diagrams in this case. In such a way, we have the 

expansion 

-RiM= hk^^ X ^ k + ^ x x L=k +... (17) 

Here we have indicated explicitly that certain index 1 < i < n refers to the outer lines 
in the n-component case. The waved line, given by (|12|) . corresponds to £ = 1. The 
combinatorial factors 4, 32, etc., are included in the diagrams not distinguishing the two 
orientations with respect to the vectors k and — k as different. The resummation of 
expansion 1)17(1 can be made by adopting our method presented by Eqs. (|13|) and (|14|). 

The Dyson equation ® can be formulated in accordance with the variational principle, 
i. e., it follows from the extremum condition dF(G)/dGi(\s) = of the reduced free energy 
function 

D(G) . (18) 



i? ( G ) = -^E ln [ 2 ^( k )] + E 



0i(k)G f (k) - \ 



Here 6>j(k) = tq — 5 • 6^1 + ck 2 and F(G) depends on the set of discrete variables Gj(k), 
as consistent with the diagrammatic definition of D{G). This function is related to the 
free energy F = —TlnZ via F(G) = F/T. First it has been obtained in |26j for the case 
n = 1. The fact that F(G) provides the diagrammatic representation of the reduced free 
energy can be proven as follows. According to @, dF(G)/dGi(k) = holds, neglecting 
the remainder term i9j(k). It means that 



?Gl(k )_V^(x» = ^ (19) 



is true within the diagrammatic perturbation theory. Besides, it is straightforward to 
check that F(G) = F/T holds at ro — ► +00. By integration in (J19|) over r$ from +00 to 
any finite value we find that F{G) = F/T is valid in general. 

4 The correlation function and susceptibility below T c 

Some important relations between the correlation function, the long-range order parame- 
ter M (e. g., magnetization or polarization), and susceptibility Xj following directly from 
the first principles, are considered in this section. 

We have defined the correlation function in the coordinate representation as 

G;(x) = ( w (xi) w (xi + x)> = V- 1 Gi(k) e 4kx ■ (20) 

k 

For simplicity, first, let us consider the one-component case. In this case (omitting the 
index 1) we have [Jj 

M 2 = lim G(x) = G(0)/V at V -»■ 00 , (21) 
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or 

G(x) = M 2 + G"(x) , (22) 
where G'(x) tends to zero if x — ► oo. In the Fourier representation (|22|) reduces to 

G(k) = 5 ki0 KM 2 + G"(k) (23) 

where G"(k) is the Fourier transform of G"(x). The susceptibility, calculated directly from 
the Gibbs distribution, is 

X = lim lim ^ = lim lim / (<5(x) - (p) 2 ) dx = G'(0) . (24) 

In this limit (</j) = M holds, the latter, however, is not correct at h = when (ip) = 0. 
The considered limit for x exists and G'(0) has a finite value in our case of n = 1, since the 
correlation function G'(x) is characterized by a finite correlation length £, which ensures 
the convergence of the integral in (|24|) . 

Consider now the case ra > 1. If an external field is applied along the i—th axis with 
i = l (even if h — > 0), the longitudinal Greens function Grii(k) = Gi(k) behaves in a 
different way than the transverse one Gj_(k) = Gj(k) with j ^ 1. It is a rigorously stated 
fact that G_l(0) diver ges as M/h if h — * below T c , which is related to the divergence 
of the transverse susceptibility in this case. In analogy to (|23|) and (|24j) we have 

G||(k) = 5 kj0 ^M 2 + Gj|(k) , (25) 
X|| (fc) = dM/dh = G\ l (0) . (26) 

Our further analysis shows that the longitudinal susceptibility X||(^0 diverges at h — ► 
for 2 < d < 4 and n > 1, i. e., Gji(k) diverges at k — > 0. Note that Gy(k) = Gy(k) holds 
at k ^ 0. The long- wavelength divergences of the transverse and longitudinal correlation 
functions below T c is known in literature as the Goldstone mode singularities established 
by the Goldstone theorem |^7| El • 



5 Generalized scaling hypothesis 

According to the known jH] scaling hypothesis, the correlation function above the critical 
point, i. e. at T > T c and T — > T c , can be represented in a scaled form 

GiM-e^gM), (27) 

where £ is the correlation length, r/ is the critical exponent, and gi(z) is a scaling function. 
Since £ ~ holds, where t = (T/T c ) — 1 is the reduced temperature, the correlation 
function can be represented also as 

G?i(k) ~ r^Qi {kr v ) , (28) 

where 7 = (2 — rj)u. Since the phase transition occurs merely at a single point h = t = 
in the h—t plane, there exists a way how the scaling relations like l|27jl or Q28|) can be 
continued to the region t < passing the singular point h = t = 0. Eq. Q27|) is not valid 
at /i = and f < in the case of n > 1, since Gj, (0) and, therefore, the correlation length 
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£ diverges at h — ► [cf. Eq. (|'24|)] . The known scaling relations are recovered assuming 
that the physical picture remains similar if we approach the critical point like t—*st and 
h — > s°"/i, where s < 1 is the rescaling factor. Thus, the distance from the critical point 

1 /2 

i = (t 2 + h?/ a ^j and the polar angle 9 = arctan ^/i 1 /""/^ in the t - h}/ a plane are two 
relevant scaling arguments. According to this discussion, a suitable generalization of the 
scaling relation (|28j) is 

G<(k) ~ i-i 9i (kt-",0) , (29) 

which is true at t — ► for any given values of kt~ v and 9. Consider Gj_(0) at a small 
negative t. Taking into account that h l / u ~ £ (7r — 0) and M oc (— t)' 3 ^ hP'°(iT — 9)~@ hold 
at — > 7r, the correct result Gj_(0) ~ M/h is obtained in this limit if g±(0, 9) oc (it — 9)~ a 
holds at 7r — 9 —* and the scaling dimension is cr = /5 + 7. By generalizing the scaling 
relation M oc (-t) 13 to M oc i* 8 (at a fixed 0) we obtain also the correct behavior M oc h 1 ^ 
at t = 0, where 5 = 1 + 7 Eq. (|2l?)) makes sense for Gii(k) at k / 0. According to l(23|) 
and (|26j). the longitudinal susceptibility is X|| = Gn (0) = G||(+0), where G||(+0) denotes 
the value of the Greens function at an infinitely small, but nonzero k value. It is easy 
to check that IJ29JI reproduces the known scaling behavior of x\\ for t > both at h = 

(X|| oc t^) and at t = (x\\ oc h?~ l ). In the limit h -> Eq. (|15J) yields 

Gi(k) ~ I t T 7 5+ (jfe I t \~ u ) at i > (30) 
Gi(k) ~ I t |~ 7 (jfe I t |~ y ) at f < , (31) 

where gf(z) = gi(z,0) and g^{z) = gi(z,n). The analysis of our diagrammatic equations 
confirm the scaling relations (|30|) and (|31|) . It shows also that, in the case of the order 
parameter dimensionality n > 1, both the longitudinal and the transverse Greens functions 
diverge at k — ► when T < T c . It means that g^(z) diverges at z — > for n > 1. In 
any case we have g^{z) oc z~ 2+ri at z — > 00, which means that the correlation function 
continuously transforms to the known critical Greens function Gj(k) ~ k~ 2+v at t — > 0. 

6 Exact scaling relations and their consequences 

In distinction to Sec. 03 here we consider other kind of scaling relations which also are 
relevant to our further analysis. These are exact and rigorous relations between the cor- 
relation function and parameters of the Hamiltonian H/T = Hq/T + Hi/T, where Hq/T 
is given by © and H\/T represents the if 4 contribution (j2j) at Uk = it. 

Following the method described in Appendix B of jS], the Hamiltonian H/T is trans- 
formed to 



H/T = (^-e-^,i + P 2 ) I *i(p) I 2 (32) 

i,p 

+ Vf 1 J2 *t(Pl)*<(P2)*j(P3)*i(-Pl -P2 -P3) , 

*J.Pl.P2,P3 

where * 4 (p) = u a c- da / 2 ipi(k), P = c 2Q u~ a k, R = r c da u~ 2a , and e = 5u- 2a c da . Here 
a = 1/(4 — d) and the summation over p takes place within p < po = c 2a u~ a A in 
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accordance with the rescaled volume V\ = Vc 2da u . According to (|32|) we have 



Gi(k) = (| <pi(k) | 2 ) =c da u- 2a (\ *i(p) | 2 ) = c^- 2 ^(p,Po, J R,e,^i) , (33) 



where, at fixed d and n, gi is a function of the given arguments only. The thermodynamic 
limit exists at k 7^ 0, so that in this case we can write 



where the scaling function g~i represents the limit of Tji. If R is varied, then the model 
with Hamiltonian (|32j) undergoes the second-order phase transition at some critical value 
R = R c (po) < 0. Thus, Eq. (|M|) can be rewritten in new variables with a scaling function 
gi(p,Po,t), where t = 1 — (R/R c ) = 1 — (ro/ro c ) is the reduced temperature, tq c being the 
critical value of tq. Thus, we have 



Based on (j34j) and (|35|) . we can make some conclusions about the scaling of the critical 
region for the reduced temperature t, as well as for the wave vector k at t = and also at 
a fixed t < 0. The latter case is relevant for the long- wave limit k — > at n > 1. By the 
critical region we mean the region inside of which the correlation function is well described 
by a certain asymptotical law. According to (|34|) and (|35j) . the width of the critical region 
t cr it or k cr n, as well as the coefficients in asymptotic expansions in powers of & at a fixed 
t (t = for n > 1 or t < for n > 1) can be written in a scaled form with (or without) 
power-like prefactors and scaling functions containing single argument po = c 2a u~ a A. 

The limit u — > is important in our consideration. To extract exact critical exponents 
from our equations, it has to be ensured that, inside the critical region, the remainder 
term #j(k) in Q is much smaller than any term in the asymptotic expansion of 1/Gj(k). 
This is possible at u — > if these expansion terms and also the width of the critical region 
do not tend to zero faster than any positive power of u 0. For the scaling functions of po, 
the limit u — > is equivalent to the limit A — > 00 at d < 4. The relevant quantities (t cr it, 
k cr u, and expansion coefficients at k powers) can tend to zero exponentially at u — > only 
if the corresponding scaling functions of c 2a u~ a A do so. The latter would mean that these 
quantities are very strongly (exponentially) affected by any relatively small variation of 
the upper cutoff parameter A at large A values. It seems to be rather unphysical, since the 
long-wave behavior at a fixed t (also at t < 0) cannot be so sensitive to small variations in 
the short-range interactions. Due to the joining of the asymptotic solutions, the fact that 
the expansion coefficients in k power series at t = are not exponentially small in u means 
also that the same is true for the expansion coefficients in | t | power series at | t |— > 0. 
Thus, only the remainder term i?i(k) tends to zero faster than u s at any s > 0, provided 
that our solution represents an analytic continuation (see the end of Sec. El) from the stable 
domain ro > where our equations have originated and where this basic property of #j(k) 
follows directly from our derivations. It implies, in particular, that the solution below 
T c should coincide with (|3*T|). as consistent with the existence of continuous second-order 
phase transition. If we have a smooth solution of this kind, then the critical exponents 
can be determined by considering suitable limits [2] (u — > and k ~ u r k cr it(u), or u — » 
and t ~ u T / u t cr it{u) with r > 0) not only at T = T c and T — > T c , but also at T < T c . In 
this case the remainder term #j(k) is negligibly small 5_. 



lim lim G i (k)=c da u- 2a g i (p,p ,R) 

8^0 V— >oo 



(34) 



9i(p,Po,R) = 9i(p,Po,t) ■ 



(35) 
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7 The low— temperature solution at n = 1 



Let us now consider the solution of our equations below the critical point starting with the 
case n = 1. The symmetry breaking term with 5 is irrelevant at n = 1, therefore we set 
5 = 0. According to (|21j) . 1/G(k) vanishes at k = in the thermodynamic limit V — > oo, 
so that the equation © (taking into account (JJJ) and omitting the irrelevant correction 
i9j(k)) for scalar order-parameter field (n = 1) can be written as 

l/(2Gi(k)) = cfc 2 + i?i(k) -i?i(0) at k/0, (36) 
l/(2Gi(0)) = r Q + 2uG + R 1 (0) = 0, (37) 



where (for arbitrary n) 



G=U(x))=V- 1 J2G l (k), (38) 



ag*(G,i) 

^ (k) "" ac<(k) ■ (39) 

To simplify the notation, further we shall omit the index i = 1 in the actual case of n = 1. 

According to (|21[) , single term with k = gives a nonvanishing contribution to (|38|) at 
V — > oo, while the contribution of all other terms may be replaced by an integral, i. e., 

G = M 2 + (2vr)- d / G'(k) dk . (40) 



Similarly, terms with M 2 , M 4 , M 6 , etc. appear in due to the contributions provided 
by zero-vectors related to some of the solid lines. For instance, the zeroth-order term 
reads 

E(°)(q) =-2n^2M 2 G(q) + (2^)- d y > G / (k)G / (q-k)dk^ . (41) 

The terms with the spontaneous magnetization M appear in our equations in a natural 
way if we first consider a very large, but finite volume V, which then is tended to infinity. 
They appear as a feedback which does not allow the right hand side of (|37|) to become 
negative, by keeping it at 1/(2(7(0)) ~ 1/V, when ro goes to large enough negative values. 

The actual model at n = 1 belongs to the Ising universality class characterized by a 
finite correlation length at T < T c . It means that G"(0) has a finite value and 1/G(k) 
transforms to zero at k = by a jump. According to ()36|) . it means that i?(+0) ^ R(0) 
holds, where the value of R(k) at an infinitesimal non-zero k is denoted by i?(+0). To 
show that this is really possible, consider the contribution (denoted by R^°\k)) of the first 
diagram in (Jllj) which yields 

J2(°> (+0) = (0) + ? ^ ^ • ( 42 ) 

1 + 4n [2M 2 G'(0) + (2?r)~ d J G" 2 (q) dqj 

From (|l2|) we see that i?(+0) ^ R(0) holds, in general, if G'(0) has a finite value. Therefore 
such a selfconsistent solution is, in principle, possible. 

Consider now the solution at ro — * — oo and small u, i. e., at low temperatures. In this 
case we find a solution such that 

M 2 = -r /(2n), G(k) = -1/ (Ar ) , E(k, C) = /(C, A) , (43) 
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hold at any fixed u > and k ^ 0, if ro — ► — oo, where A is a constant independent of 
ro, c, and A, /(£, ^4) is a function of the given arguments. One expects that A tends to 
some universal constant at u — > 0. In this case, at ^4 = 4, our solution coincides with the 
Gaussian approximation 

G'(k) = — — . (44) 

v ' -Ar + 2ck 2 v ; 

The correction ~ ck 2 has been neglected in (|43|). 

Condition M 2 = —ro/(2u), corresponding to the minimum of 0, holds for any physical 
solution if ro — ► — oo since fluctuations are suppressed. This follows from (|3*7|) and (|40|). if 
the main terms ro and 2tt M 2 are retained in ()37jl . 

Quantity S(0, C) i n the denominator of (|12|) diverges in the thermodynamic limit 
V — > oo like V or even faster, i. e., like with \i > 1, as shown in Appendix. The 
divergence appears due to the special contributions provided by zero wave vectors k = 
related to some of the solid lines in the diagram expansion Ijlljl of S(0, Q at T < T c . 
These diverging terms make the analysis below T c more difficult as compared to the case 
T > T c discussed in 0. 

In spite of the divergence of S(0, £), a single term with q = in (jSJ does not contribute 
to D*(G, C) and R(k) if V -> oo. Really, if E(0, () diverges as - (V/V ) s (V Q is the volume 
of an elementary cell) with s > 0, then this single term is approximately — (s/2) -In (V/Vq), 
whereas the whole sum is proportional to V/Vq. On the other hand, S(0, Q appears in the 
denominator of the corresponding term if the derivative with respect to G(k) is calculated 
in (JHJ), therefore, this term cannot be by a factor V larger than other terms, i. e., it cannot 
give an especial contribution. 

The perturbation sum for S(k, Q at k ^ also contains terms diverging at V — > oo. 
In a normal case, the constraint k = for wave vectors of m solid lines in a diagram means 
removal of m integrations over wave vectors. However, for some distributions of the zero- 
vectors this condition is violated, which yields the diverging terms. We have shown in 

o y/A 

Appendix that the divergent terms contain insertions like, e. g., n V/X n with 2m 
outer solid lines having fixed k = vectors. As shown in Appendix, a resummation of 
these insertions gives a non-divergent result. Besides, the simplified analysis which ignores 
these insertions is qualitatively correct, as regards the general scaling form of the solution. 
However, for a complete formal correctness we should take into account the fact that 
specific values of scaling functions can be renormalized by these zero-vector-cumulant 
(see Appendix for explanation) insertions. Further we shall call the terms without such 
insertions the "normal" terms or contributions. Our results at T < T c are completely 
consistent with those at T = T c and T > T c and agree with the non-perturbative RG 
analysis provided in Sec. 19.51 It confirms the statements made in this paragraph. 

An important property of the "normal" contributions to S(k, C) at k ^ is that any 
term, where zero wave vector is related to a waved line, vanishes in the thermodynamic 
limit V — > oo. It holds because the waved line vanishes due to the divergence of S(0, £). 

Condition E(k, £) = /(£, A) at k ^ holds if ro — > — oo because E^(q) and all terms 
of the sum (|14|) in this limit depend merely on parameters £ and A for any fixed u. The 
latter is true in the approximation where no zero-vector-cumulants are included, since 
the main terms come from the diagrams in ()11|) if we extract the contribution, containing 
no integrals, where one half of the solid lines have zero wave vectors k = (it yields a 
factor M 2 for the first diagram, M 4 - for the second diagram, and so on). Besides, the 
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waved lines should have nonzero wave vectors, as explained before. According to (|43|) each 
replacement of G(0)/V = M 2 by an integral over wave vectors produces a factor ~ r^ 2 . 
This is the reason why the main terms at rg — ► — oo contain no integrals. 

These main terms lead to asymptotic solutions (at ro — * — oo) for quantities 
<9£(q, C,)/dG(h) and dD*(G, ()/dG(k) at q = k 7^ in the form where r is multiplied by 
some function of £ and A, namely, i?(k) = — Q(A) tq. This main contribution vanishes at 
k = 0, since the derivation dT,(0, £)/dG(0), implicated in the calculation of R(0), means 
removal of a diverging factor G(0) ~ V. Quantity S(0, Q diverges, whereas the waved 
lines in this case contain vanishing factors ~ 1/S(0, £) providing finite value of the deriva- 
tive at V — > 00. Therefore, by substituting i?(k) = — Q(A) ro into (|36|) and retaining only 
the leading terms, we obtain a selfconsistent equation 2Q(A) = A for the unknown ampli- 
tude A. Following the consideration in Appendix, the inclusion of zero-vector-cumulant 
insertions can merely renormalize the function Q(A). 



8 Asymptotic solution at T — > T c 

Our equations provide the asymptotic solution at T — > T c , where T < T c , in the scaled 
form (|31[) which allows a unified description provided here for n > 1. 

Like in [5J, here we assume that vq is the only parameter in © which depends on 
temperature T and the dependence is linear. At the critical point ro = ro c we have 
l/Gj(+0) = for all i at lim lim , so that the Dyson equation (J0J in this limit reads 

5^0 V^oo 



1 dr 



2Gi(+0) dT 
1 1 



• A + 2« (G-G*) +Ri(+0)-R*(+0) (45) 
ck 2 + Ri(k) - Ri(+Q) . (46) 



2Gi(k) 2Gi(+0) 

Here A =| T-T c \ is an analog of | t \ in (|3T|) . G* and R*(+0) = R*(0) are the values of G 
and Ri(0) calculated at the critical Greens function Gj(k) = G*(k) considered as a known 
fixed quantity. Due to the symmetry breaking term 5 ■ 5n in (j3J), only the longitudinal 
component 1/Gi(0) = 1/G||(0) becomes as small as ~ 1/V when ro is decreased below 
roc, giving rise to the magnetization M 2 = G\\(0)/V. According to (|46l) . the condition 
1/G||(0) = at V -> 00 means 

1/ (2G|| (+0)) = R\\ (+0) - R {{ (0) . (47) 

Equations Q45JI and are analogous to (48) and (49) in 5, derived for T > T c , 
and similar method of analysis is valid. Namely, correct results for the Greens function 
within the asymptotical region k ~ l/£ can be obtained considering the limit u — > 
and A ~ u r ^ u A cr i t (u), and formally cutting the integration over Gj(k) and G*(k) in (|45|1 
and by fc < A = ii -r £;~^(u)/£. Here A cr ^(u) and k cr i t (u) are the widths of the 
critical regions for A and k, respectively, r is any positive constant, and £ is an analog of 
the correlation length. According to the generalized scaling hypothesis in Sec. [51 one may 
set £ = £(T = T + A). Note that our equations (f4"5|) to (jlTjl do not contain ro c . Only 
such a form is acceptable in this analysis: contrary to the critical exponents the critical 
temperature is essentially affected by the short-wave fluctuations. 
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We start our analysis with the "normal" (not diverging at V — > oo) terms discussed in 
Sec. [7] by using the well known scaling relation 



2/3 = (d - 2 + t})v = dv - 7 



(48) 



relevant to M ~ A@. Our diagram equations become exact in the theoretical limits 
considered, and it allows in principle to find the exact critical exponents by solving these 
equations (SecEJ). Therefore, the asymptotic solution satisfies the existing exact relations 
between the critical exponents, including (48). In fact, (|48j) is the necessary condition at 
which the diagrams below T c have certain scaling properties, similar to those above T c , 
where the exponents are common for the diagrams of all orders. 

Selfconsistent solutions for E(q, C) and for dD*(G, ()/dGi(k) can be found at 2 < d < 
4, having similar form as in the case T > T c 0. This is true for the "normal" terms 
because partial contributions in p^p proportional to M°, M 2 , M 4 , etc., (corresponding 
to cases where zero-vectors are related to 0, 1, 2, etc., solid lines in diagrams pip) have 
the same form and the same common factor /±~ 2 "t+ du as at T > T c . It is easy to prove 
this property for diagrams of all orders by a simple normalization of all wave vectors to 
A u , like q' = q/A u , with the cutting of integration discussed above. Contrary to the case 
T > T c , the main term of <9E(q, ()/dGi(k) at T < T c for given u, d, and n has the form 



obtained by the above normalization, where k' = kA ~ v . The additional (second) term 
is due to the following. If in some diagrams of pip the wave vectors for some solid lines 
are fixed ki = (this yields an especial contribution merely at T < T c ), then the wave 
vectors of some other solid lines also have a fixed value k2 = q. This produces the second 
(extraordinary) term after derivation of the lines with a fixed wave vector q. The simplest 
example is provided by the first diagram in pip represented by (|4ip. The first term 
in (|4ip . which at arbitrary n reads — 4uM 2 Gi(q), yields an extraordinary contribution 
— 4u5j i i5 q i c M 2 to l|49p. In general, any diagram of pip can be represented as two separate 



gives an extraordinary contribution provided that no more than one pair of the connecting 
lines meet at the same node like r ^x<^ . Since solid lines always make closed loops, m 
is even number. The same index i is related to all lines of one loop, therefore, even 
number of lines among those connecting the blocks can be associated with i = 1. Thus, 
an extraordinary contribution is provided by fixing i = 1 and ki = for m — 1 connecting 
lines, since this implies the constraint i = 1 and k2 = q for the remaining one line. 
However, only those configurations give a nonvanishing contribution at V — > oo where no 
more than one wave vector ki = refers to each of the configurations r ^-<^ . In the 
opposite case the waved line also has zero wave vector yielding vanishing contribution due 
to the divergence of S(0, (). 

The additional term in (|49p essentially differs from the ordinary (first) term by factor 
5^qVA du because in this case the derivation procedure does not mean removal of integra- 
tion over wave vectors. Both terms provide dD*(G, Q/dGiQs.) = A 7 r(k',£), single term 
with q = in jSJ) being negligible, as discussed in Sec|7| This scaled form is valid also for 
i > 1 and agree with that at T > T c derived in The term with <5k, q does not give a 
contribution to dD*(G, Q/dGi(k) at V — * oo if k = 0: any of the partial contributions to 
S(0, £), where zero-vectors are related to the connecting solid lines in a configuration like 



0E(q, C)/dGi(k) = V- l A^ Yiiq', k', () + <V4, q A dl ^ F(q', k' 







(49) 



parts connected by m solid lines, e. g., like 
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, is by a factor ~ Gi(0) = Gm(0) ~ M 2 V larger than the corresponding 



term of 9S(0, Q/dGi(0), where factor G||(0) is removed according to derivation of the 
connecting lines with k = 0. No diverging factor is removed at k = +0 if G||(+0) has 
a finite value, therefore, i?(+0) and R(0) are not identical at n = 1. In any case, this 
peculiarity refers only to the longitudinal component, so that R±(+0) = R±(0) holds in 
general. 

The above discussed peculiarities with the extraordinary terms are irrelevant, as re- 
gards the exponents in the asymptotic expansion of Gj(k). The magnetization M has 
corrections to scaling with the same exponents as other terms in 1)45 [I and (|46|) (although 
some expansion coefficients can be zero). This is consistent with the additional condi- 
tion (|47|) - the equation for M. The corrections to scaling have the same origin at T < T c 
and T > T c : the term ck 2 in (|46|) is by a factor ~ A 2 ^ -7 smaller than the leading term and 
the term "1" in (|12|) gives (in the same sense) a correction ~ A 27_dl/ . In such a way, the 
general scaling form of the solution for Ri(k) and Gj(k) appears to be the same at T < T c 
and T > T c . Following the Appendix, it remains unchanged in the complete analysis which 
takes into account all possible distributions of k = vectors in the diagrams of (jll|) . and 
not only the "normal" contributions without the zero-vector-cumulant insertions. 

Thus, for the spational dimensionality 2 < d < 4, the correlation function at T < T c 
has similar singular structure as in the case of T > T c 0, i. e., we have an asymptotic 
expansion 

G J (k) = ^A^ 9 f ) (kA^) (50) 

with 70 = and correction exponents 

7£ = n e {2u - 7) + mi (27 - dv) (51) 

valid for £ > 1, where rig and are integer numbers > 0, and + > 1. Therefore, 
using the same arguments as in the case of T > T c [Sj, we conclude that the possible values 
of exponents 7 and v are given by ()15|) and 1)16(1 . It agree with the generalized scaling 
hypothesis in Sec. El which tells us that the values of the exponents must be the same at 
T>T C and T < T c . 

9 The asymptotic long— wave behavior below T c at n > 1 
9.1 Discussion of the existing results 

According to the conventional believ jT^l [TTl [THl H3 HOI HH H31 HU HS1 12H] , the transverse 
Greens function Gj_(k) diverges like k~ x± with Aj_ = 2 at k — > below T c for the systems 
with 0(n > 2) rotational symmetry. It corresponds to the Gy(k) ~ k d ~ 4 divergence of the 
longitudinal Greens function. Besides, the singular structure of the correlation functions 
is represented by an expansion in powers of k 4 ~ d and k d ~ 2 P21l23j. Formally, our results 
agree with these ones at Aj_ = 2. Nevertheless, below we will show that Aj_ < 2 holds near 
two dimensions at n = 2. 

As usually accepted in lattice models, here (in this subsection) we define that all the 
parameters of the normalized Hamiltonian H/T © are proportional to the inverse tem- 
perature 1/T. In this case tq is negative. We propose the following argument. The 
assumption that Gj_(k) ~ a{T)k~ 2 holds (with some temperature-dependent amplitude 
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o(T)) in the stable region below the critical point, i. e., at T < T c /C, where C is an arbi- 
trarily large constant, leads to a conclusion that the critical temperature T c continuously 
tends to zero at d — > 2 (supposed d > 2). Really, at A^ = 2 we have 

^(x)) = ^£G,(k)^M 2 (52) 



-(2„)- 



where S'(d) is the area of unit sphere in d dimensions. Since the amplitude of the transverse 
fluctuations never can vanish at a finite temperature, Eq. 1)52(1 implies that the average 
(</? 2 (x)) diverges at T = T c /C when d — > 2, if T c remains finite. Thus, we obtain an 
unphysical result unless the critical temperature T c and, therefore, a(T c /C) tend to zero 
at d -> 2. 

On the other hand, it is a rigorously stated fact |29U3U| that the classical 2D XY model 
undergoes the Kosterlitz-Thouless phase transition at a finite temperature T KT . It means 
that a certain structural order without the spontaneous magnetization exists within the 
temperature region T < T KT . There is a general tendency of disordering with decreasing 
the spatial dimensionality d, and not vice versa. Thus, since the structural order exists at 
T < T KT and d = 2, some kind of order necessarily exists also at T < T KT and d > 2. Since 
the classical XY model undergoes the disorder — > long-range order phase transition at 
d > 2, this obviously is the long-range order. Thus, the critical temperature at d = 2 + e 
is T c > T KT / for an infinitesimal and positive e, and it drops to zero by a jump at 
d = 2 — e, as consistent with the rigorous consideration in [2D] . 

The classical XY model belongs to the same universality class as the actual if 4 model 
at n = 2, which means that both models become fully equivalent after a renormalization 
(a suitable renormalization will be discussed in Sec. 19. 5j) . Thus, T c does not vanish at 
d — * 2 (for d > 2) also in the (p 4 model. In such a way, the assumption Gj_(k) ~ cl(T) k~ 2 
leads to a contradiction. In the stable region T < T c /C, the Gaussian approximation 
Gj_(k) — l/(2cfc 2 ) makes sense at finite not too small values of k. The above contradiction 
means that the Gaussian approximation with Aj_ = 2 cannot be extended to k — * in 
vicinity of d = 2. The contradiction is removed only if Aj_ < 2 holds at d — >• 2 in the actual 
case of n = 2. 

It has been stated in [231 I25| that the essentially Gaussian result Aj_ = 2 of the 
perturbative RG theory should be exact. However, some of the obtained "exact" results 
are rather unphysical. In particular, we find from Eq. (3.6) in [21] and from the formula 
(it 2 ) = — 6A/uo given in the line just below that 

<^ 2 (x)> = (iV-l) f (53) 

holds, where vr(x) is the transverse (N — l)-component field which, in our notation, is 
composed of n— 1 components labeled by index j ^ 1. Eq. (|5^|) represents a senseless result, 
since (7r 2 (x)) given by this equation diverges at d — * 2. It is clear that (7r 2 (x)) cannot 
diverge in reality, as it follows from the Hamiltonian density (2.1) in [21] (Hamiltonian Jy) 
in our paper): any field configuration with diverging vr 2 (x) provides a divergent 7r(x)- 
dependent term 

~ \ | W(x) | 2 +|V(x)) 2 (54) 
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in the Hamiltonian density and, therefore, gives no essential contribution to the statis- 
tical averages. The result f)53f) corresponds to a poor approximation where the second 
term in (|54|) is neglected. In a surprising way, based on Ward identities, authors of [21] 
and related papers have lost all the purely transverse diagrams, generated by the term 
(7r 2 (x)) 2 , and stated that this is the exact result at ro — ► — oo as well as at k — > 0. 
According to |23l I24| , the actual transverse term appears to be hidden in a shifted lon- 
gitudinal field s, which is considered as an independent Gaussian variable (cf. Eqs. (3.5) 
and (3.6) in 23 J. Obviously, this is the fatal trivial error which leads to the above dis- 
cussed unphysical result (|53j) . since the determinant of the transformation Jacobian (from 
7r, s to 7T 3 s) is omitted in the relevant functional integrals. In this manner the ip 4 model 
can be immediately reduced to the ip 2 model by considering s = <p> 2 as a new variable! 
Including the Jacobian of the nonlinear transformation, the resulting model, however, 
remains non-Gaussian. Since (|53|) comes from 

(vr 2 (x)) = (N - 1) (2iry d J G±(k) d d k , (55) 

the unphysical divergence of (7r 2 (x)) means that the predicted Gaussian form of the trans- 
verse correlation function G±(k) is incorrect. Another aspect is that the method used 
in [2H1 I2U gives X± = 2 also at n = 2 in contradiction with our previous discussion 
concerning the known rigorous results for the XY model. 

Our consideration does not contradict the conventional statement (see [23 and refer- 
ences therein) that the Gaussian spin-wave theory becomes exact at low temperatures, 
but only in the sense that it holds for any given nonzero k at T — > 0, and in the limit 
lim lim in particular. However, the actual limit of interest is k — > or, equally, lim lim. 

Therefore, contrary to the assertions in |14[ I25j . it is impossible to make any rigorous con- 
clusion regarding the exponent Aj_ (or any related exponent) based on the fact that the 
Gaussian spin-wave theory becomes exact at T —* 0. One has to prove that the limits can 
be exchanged! There is a reason to believe that the limits cannot be exchanged first of all 
because the critized here treatments with exchanged limits lead to unresolvable problems 
at d — > 2. There is also a well studied example - the XXX quantum spin chain, where it 
is straightforward to see that the distance — * oo and T — ► limits cannot be exchanged 
in the correlation function considered there |31j . 

This problem persists in the classical treatment of the many-particle systems JH] 
which, in essence, is based on the Gaussian spin wave theory at T — > 0, as well as in the 
hydrodynamical description in |14j . where the known results of the Gaussian spin wave 
theory [Eqs. (5.1a) and (5.1b)] have been implemented for a complete description. The 
treatment of ^S], evidently, is not exact, since it breaks down at d — * 2 for the two- 
component (n = 2) vector model (where T c remains finite and we fix the temperature 
< T < T c ) just like we have discussed already - the average (7r 2 (x)) is given by the 
integral (6.8) in j£] which diverges in this case (supposed (2ir)- 3 d 3 k -> (2ir)- d d d k). 

A slightly different perturbative RG approach has been developed in [2J to analyze 
the nonlinear a model. In this case the modulus of <p(x) is fixed which automatically 
removes the divergence of (7r 2 (x)). A finite external field h has been introduced there to 
make an expansion. The correlation functions have the power-like singularities of interest 
only at h = +0, which means that in this case we have to consider the limit lim lim, i. e., 

the limit h — > must be taken first at a fixed nonzero k (p in formulae used in (U). The 
results in [2J are n °t rigorous since the expansions used there are purely formal, i. e., 
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they break down in this limit. Besides, contrary to the approximations in |21j . it should 
be clear that the exact renormalization is a rather nontrivial problem which cannot be 
reduced to a finding of only two renormalization constants. If, e. g., we make a real-space 
renormalization of the Heisenberg model, say, with the scaling factor s = 2, then the 
statistically averaged block-spins of the Kadanoff transformation (composed of s d original 
spins) do not have a fixed modulus - simply the original model does not include the 
constraint | ty?(x) |= const for the block averages. It means that the transformation with 
any finite s yields a Hamiltonian of form different from the original one, i. e., the original 
Hamiltonian with merely renormalized coupling constant can never be the fixed-point 
Hamiltonian. 

Another approach, which is based on effective Lagrangians, has been developed in |28j . 
However, due to several rough (fatal) errors in Sec. 3 of this paper we cannot appreciate 
the basic results. For example, we have found that the transformation 

v% = riv^tl' 1 , H n = Off (56) 

brings the Lagrangian 

C = ±DjDj+ \m 2 U+ \\{U? - H$ (57) 

with 

ZV(x) = fy£(x) + <v(x)0(x) (58) 
back to its original form after a space-independent rotation f2. Namely, 

c(y%,H n ,n$) =c(v^Hj) (59) 

holds, as it can be verified by a simple substitution taking into account that Q,A VLB = AB 
holds for any vectors A and B, since the rotation of the coordinate system does not change 
the scalar product, as well as that d^Qcj) = ^ d^cj) obviously is true for a space-independent 
matrix f2. As a consequence, the partition functions obey the equation 

z(v%,H n )=z(v IM ,H) (60) 

with v% and H n defined in But Eq. (3.7) in [2S| disa grees with (|56|): it contains an 
odd term Qd^Q^ 1 which reduces to in the case of space-independent Q. The following 
relevant equations on page 247 contain similar errors, i. e., it is easy to verify that they do 
not hold in the simplest case of a space-independent rotation. Some numerical support 
of this theory can be found in literature |321 I33j . where a finite-size scaling within a 
transient region of very small fields comparable with L -3 (where L is the linear size of 
3D lattice) have been considered. Our theoretical predictions, however, refer to the limit 
lim/j^o nm L^oo ) therefore the tests made in |32l I33j are of little interest here. 

Concluding this subsection, it is worthy to mention that the experimental measure- 
ments of susceptibility \ depending on field h in isotropous ferromagnets like high-purity 
polycrystalline Ni [31] are incompatible with the conventional RG prediction x ~ 
Moreover, our Monte Carlo results for 3D XY model discussed in Sec. are incompati- 
ble with this prediction, as well. Thus, according to the given theoretical arguments, the 
conventional claims that the Gaussian approximation is asymptotically exact at k — > 
simply cannot be correct, and there exist also numerical evidences for this. 
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9.2 The leading asymptotic behavior 

Let us now discuss the solution below T c at small wave vectors k —* within our diagram- 
matic approach. By analyzing several possibilities we have arrived to a conclusion that 
the true physical asymptotic solution for 2 < d < 4 and n > 1 is 

Gj_(k) ~ a*T A ^ , Gy(k) ~ 6AT A ll (61) 

whith some amplitudes a and b, and exponents 

d/2 < A ± < 2 and Ay = 2A_l - d . (62) 

Only in this case the exponents on the left hand side of equation (JBJ) for 1/Gj(k) coincide 
with those on the right hand side, if calculated by the method developed in for the 
case T = T c . Besides, only at d/2 < Aj_ < 2 we arrive to a solution which coincides both 
with the scaling hypothesis 1(31)1 and general renormalization group arguments discussed 
further in Sec. 19.51 

Below we show that (|61|) and (|62jl really represent a selfconsistent solution of our 
equations. According to (|62|) . An > holds, so that Gy(k) diverges at k — > 0. Thus, we 
have l/G||(+0) = 0. Eq. (@7J) then implies that R\\(+0) = R\\(0). The latter relation 
always is true for the transverse components, as already discussed in Sec. El Our analysis 
is based on Eq. I|46|). which at these conditions reduces to 

l/(2Gi(k)) = ck 2 + Ri(k) - Ri(0) . (63) 

First let us consider only the "normal" contributions without the zero-vector-cumulant 
insertions, as explained in Sec. and Appendix. In this case (at the condition (|62jl). 
the main terms in the asymptotic expansion of S(q, £) (at q / for any given A/q 
considered as independent variable) are represented by partial contributions, coming from 
all diagrams in ((lip, where either all Nj solid lines of the j-th diagram are associated 
with the transverse components Gj_(k), or m of 2m < Nj solid lines which are associated 
with Gy(k) have zero wave vector. It is true at 2 < d < 4 and (|62|) . since other terms 
provide only small corrections, as discussed further in Sec. I9..3I Since 1/S(0, () = holds 
at V = oo, only those configurations give a nonvanishing contribution where nonzero 
wave vectors are related to the waved lines. Therefore also maximally one half of all lines 
associated with the longitudinal component i = 1 can have zero wave vector, as regards 
the nonvanishing (at V — * oo) terms related to S(q, £) with q / 0. This holds because 
solid lines make closed loops and maximally each second line of a loop can have zero wave 
vector, provided that all waved and dashed lines, connected to this loop, have nonzero 
wave vectors. 

It is suitable to represent the amplitude b in (|61|) as b = b'-a 2 /M 2 . Then, by normalizing 
all wave vectors to the current value of q, we find that the selfconsistent solution of (|13|) 
has the scaled form 

S(q, C) = a 2 q d - 2X ^ (A/q, C, X±,b', u, d, n) (64) 

which is proven by the method described in detail in [S], taking into account all diagrams 
of ((11(1 and only the main terms of the asymptotic expansion at k — > 0. These are the 
partial contributions discussed in the paragraph above. It is supposed also that term " 1" 
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in (|12[) . providing a small corection at An > 0, is neglected. In the same manner, by 
normalizing all wave vectors to the current value of k, we arrive to the scaled form 

8Z(ci,()/dGj(k) = V~ 1 ak- x± Y ± (c l /k,A/k,(,X ± ,b',u,d,n) : (65) 
0E(q,C)/0Gi(k) = V~ 1 M 2 k~ d Y\\ (q/fc, A/k, (, Xj_,b', u, d, n) 

+ 8< l!h M 2 Y ll (q/k,A/k,(,\ ± ,b , ,u,d,n) , (66) 

where the term with Kronecker's symbol appears due to the extraordinary contributions 
discussed in Sec. [SJ By virtue of (|66|). and (jSJ), l|3*9")) we obtain also 

i?_i_(k) = a' 1 ^ ^ (A/k,X ± ,b',u,d,n) (67) 
i?ll(k) = ZT 1 ^!! 0|| (A/A;, A ± , u, d, n) . 



The fact that i?j_(0) is constant means that there exists the limit lining R±(k) = R±(0), 
which, according to (joTj) . implies that aA~ x± R±(0) does not depend on a and A and thus 
the scaling function (ft can be represented as 

4> ± (A/k, Aj_, b' , u, d, n) = aA' x± R_i_(0) ■ (A/k) x± +1/ [2<p ± (A/k, X ± ,b', u, d, n)] (69) 

with another scaling function <pj_ instead of (f)±. Analogous equation for 0ii reads 

011 (A/k, X±,b',u,d, n) = 6A" A iii?||(0) • (A/k) x w +1/ (A/k, X ± , b', u, d, n)l . (70) 

Note that we always consider Aj_, but not An, as an independent argument, therefore some 
asymmetry appears in formulae. By substituting Eqs. (ffiTj) to (f70|) and (|61j) into (j63|) . and 
neglecting the correction term ck 2 , we obtain 

G ± (k) = ak- x± = a(j) ± (A/k,X ± ,b',u,d,n) k~ x± (71) 
G||(k) = bk- x w = b<t>\\ (A/k,X ± ,b',u,d,n) k~ x u . (72) 

Up to now we have neglected the zero-vector-cumulant insertions in the derivation of 
Eqs. (|71|) and (|72|) . However, following the method in Appendix, these insertions can only 
renormalize the scaling functions 4>± and cj>u, so that Eqs. (|71|) and 1)72(1 hold. 

As discussed in Sec. H3 the limit u — > with simultaneous tending of to zero like 
fc ~ u r k cr it(u) (where r > 0) has to be considered to ensure correct critical exponents. 
The existence of the solution for 1)71)1 and (|72|) implies the existence of the corresponding 
limits for (j)± and (f)\\ . Note that these functions do not contain n-dependent factors in our 
scaling analysis (which is true also for <j> in Eq. (42) of jSJ), and only weak ti-dependence 
(at u — > 0) can be induced by the integration limits in (|1H[>. Thus, Eqs. ()71() and ()72*|) yield 

lim <j> ± (Au' r k-\ t (u), X x ,b', u, d, n) = B ± (X ± , b',d,n) = l (73) 

u— >0 V / ' 

lim 0|| (Au- r k-\ t (u),X ± , b' , u,d,n) = B\\ (A ± , b', d, n) = 1 . (74) 

Eqs. 1)73)1 and (|74)1 can be, in principle, solved with respect to Aj_ and b' . It follows herefrom 
that not only the exponent Aj_, but also the quantity b' = bM 2 /a 2 is universal, i. e., 
dependent only on the spatial dimensionality d and dimensionality of the order parameter 
n. In general, no universality of amplitudes is expected, so that the latter rather surprising 
conclusion refers only to the actual limit u — > 0. Nevertheless, the universality of bM 2 /a 2 
coincides with some general non-perturbative renormalization group arguments discussed 
in Sec. 19.51 which show that the restriction to small u values is purely formal. 
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9.3 Corrections to scaling 



In Sec. 19,21 we have considered only the dominant terms in the asymptotic solution at 
k — > 0. Now we shall discuss corrections to scaling. There are following sources of 
corrections. 

(i) Since R±(k) - Rj_(0) oc a~ l k Xx and A_l < 2 hold, the term ck 2 in Eq. (fo^l 
for the transverse components i / 1 produces a correction which is by factor 
£\ (k) oc a k 2 ~ x± smaller than the main term at k — > 0. In analogy, a correction 
s'i(k) oc bk 2 ~ x W is generated in the same equation for i = 1. 

(ii) According to (|64j) and 1)62 jl . at any given A/q, the term "1" in Eq. ()12j) represents a 
small correction oc a _2 g A H to the amplitude of the main term. Finally, it generates 
an amplitude correction £2(k) oc a~ 2 k x ^ in the asymptotic expansion of G(k). 

(iii) Consider partial contributions to S(q, £), coming from all diagrams in (jllj) . where 
less than one half of the solid lines which belong to loops with i = 1 have zero 
wave vector. As compared to the dominant contributions discussed in Sec. 19.21 they 
generate small corrections represented by an expansion in terms of £ 3 (fc), where 
£3 (k) oc (b/a) /c Aj -~ A H corresponds to a replacement of G±(k) with G|j(k) for one 
solid line with nonzero wave vector. 

Note that the corrections £e(k) are small at k — » only for d/2 < Aj_ < 2, so that our 
analysis cannot be formally extended outside of this interval. We have included an explicit 
dependence of £i(k) on the amplitudes a and b, since their singular behavior is relevant for 
joining of the asymptotic solutions at T — > T c . Since £[(k) oc £i(k) 63(h) holds, we have 
no more than three independent correction sources. 

The expansion in powers of £e(k) is acompanied by scaling functions depending on A/h. 
Like in (|73jl and (|74j). these scaling functions can be replaced by amplitudes which are 
independent of k, when considering the limit u — > and k ~ u r k cr it(u). This replacement 
is analogous to that at T = T c and has the same motivation [5] . It results in the asymptotic 
expansion for the Greens function 

Gi(k) = J2h(i)k- x ^\ (75) 

where relations Aj(0) = Am, 6j(0) = b hold for i = 1, and Aj(0) = Aj_, 6j(0) = a - for i ^ 1. 
A term with £ > 1 represents a correction of order e™ 1 ^ £^ 2 ^ £^ i<y ^ with the exponent 

Ai(*) = Ai(0) - m(£) ■ (2 - A ± ) - n 2 (£) ■ Ay - n 3 (£) • (A ± - Ay) , (76) 

where nj (£) > are integer numbers such that rij (£) > 1 . Note that we always allow 
a possibility that some of the expansion coefficients, in this case some of b{(£), are zero. 
The expansion in powers of £ 2 (k) oc k A ~ d and £3 (k) oc k d ~ 2 , proposed by the perturbative 
RG theory [UlESl, is recovered by formally setting A^ = 2. 

9.4 Joining of asymptotic solutions 

Consider now how our expansion (|75|) coincides with (jUlj) and (|5()j) when approaching the 
critical point. Retaining only the leading terms, the consistency is ensured if the scaling 
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functions of the dominant terms behave as g±(z) oc z x± (for i ^ 1) and g\\(z) oc z A|1 
(for i = 1) at z — ► 0, and a oc £ A ~ A± , 6 oc £ A ~ A H hold for the amplitudes in (|61|) at A — » 0, 
where A = 7/V = 2 — 77. Here £ = A - " is an analog of the correlation length. This is 
a property of the solution exceptionally at d/2 < X± < 2 that the consideration of the 
long-wave limit does not provide any constraint for the amplitude a in (|71[) . whereas the 
other amplitude b is related to a via b = b' • a 2 /M 2 . Taking into account the scaling 
law (|48|). the relations a oc £ A ~ A ^ and 6 oc £ A ~ A H mean that bM 2 /a 2 is constant at A — > 0. 
It is consistent with the statement in Sec. 19.21 that 6' = const at u — ► for any T < T c . 

The expansion (foT)|) with the exponents 7^ (|5Tj) completely agree with ((75)) and ifTSl) 
provided that the scaling functions have an asymptotic expansion 

9f\z)=Y, B t ) *- Xi{i) ( 77 ) 

at z — » 0. In this case the amplitudes ~ A-7+7«+^W have corrections to scaling 
where the main term is multiplied by oc A 7m . 

9.5 Non— perturbative renormalization group arguments 

The relation A11 = 2Aj_ — d as well as the universality of the ratio bM 2 /a 2 have a simple 
interpretation in view of some renormalization group (RG) analysis. Our ip 4 model can be 
formulated on a discrete lattice, representing the gradient term by finite differences. At 
h = +0, we consider the transformation consisting of 

(i) Kadanoff transformation replacing single lattice spins by block-spins, where each 
block-spin is an average over s d spins; 

(ii) shrinkage of the new lattice s times to return to the initial lattice constant. In 
distinction to the standard renormalization, we do not rescale the field <p. 

The distribution over block-spins is described by new Hamiltonian T S H, where the no- 
tation T s is used to distinguish from the standard RG transformation R s . The Kadanoff 
transformation does not change neither the magnetization nor the long-distance behav- 
ior of the real-space Greens functions Gj_(x) = (</?j(xi)</>j(xi + x)) = dx A± ~ d (i ^ 1) 
and G||(x) = (yi(xi)yi(xi +x)) — M 2 = bx x ^~ d at x — > 00, where a = c a ■ a and 
b = ci,-b are the amplitudes. The proportionality coefficients c a and q, are universal, since 
Gj_(k) — ak~ x± and Gy(k) ~ bk~ x ^ are the Fourier transforms of Gj_(x) and Gy(x), re- 
spectively. Taking into account the shrinkage of the lattice at step (ii), magnetization M 
is invariant of the transformation T s , whereas the amplitudes rescale as a(s) = a(l) • s x± ~ d 
and b(s) = b(l)-s x ^~ d . The modulus conservation principle is true at large renormalization 
scales s, since the variation of average modulus for large blocks of spins is related to a much 
greater increase in the systems energy as compared to a gradual long-wave perturbation 
of spin orientation. The validity of this principle is restricted by the condition that the 
mean amplitude of the relative fluctuation of the modulus has to be much smaller than 
the mean squared fluctuation of the orientation angle for the block-spins of the Kadanoff 
transformation. The non-perturbative renormalization group arguments given below are 
in agreement with our foregoing diagrammatic analysis, assuming that this condition is 
fulfilled for large s in the actual case of 2 < d < 4. Thus, the renormalized Hamiltonian 
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can be written as 

T s (H/T) ~ ^ ' 



X 



<*(l ¥>( x ) I -¥>o) 



+ Q{a- 1 / 2 ( S ) (/3± (x)} 



(78) 



at s — > oo, where the first term represents the modulus conservation principle allowing only 
those configurations with non-diverging Hamiltonian where | <^( x ) |= (fo, and Q is some 
functional of the configuration of the transverse order parameter (n — 1 component vector) 
field y?j_(x). In this case only the infinitely small (at s — > oo) transverse components <Pi(x) 
with i 1 are independent variables, since 

n 

<^o-^i(x)^5> 2 (x)/(2M) ( 79 ) 

i=2 

holds according to | (f(x) | = Not loosing the generality, we have considered the 
spatial configuration of the normalized transverse components a _1 / 2 (s) </?j_(x) as an ar- 
gument of the functional Q. According to the definition of Gj_(x), the function /j_(x) = 
(|a _1 / 2 (s) <£j(xi) • a _1 / 2 (s) <^i(xi + x)^> = d~ 1 (s)G±(~x.) with i / 1 has a universal asymp- 
totic behavior /j_(x) = x x± ~ d at x — > oo. Since this average is composed of arguments 
of Q, the sufficient condition for its universal asymptotic behavior is the universality of 
the functional Q{z(x)}. The latter is consistent with the idea that the transformation 
of Q (assuming that at any s the transformed Hamiltonian can be approximated by Q78|) 
according to some a priori defined criterion) has a fixed point 

Q*{z(x)} = lim T s Q{z(x)} (80) 

s — >oo 

which, however, may be different for each universality class. In the conventional RG trans- 
formation R s the field would be rescaled as p± 

( x ) s (d-A±)/2 ^ ^(x), so that Q in (|7H|) 
would contain no explicit scaling factor s. Nevertheless, we prefer our notation, since it 
is suited to express the modulus conservation principle. Accepting (JSUJ), any statistical 
average composed of arguments a _1//2 (s) (^j_(x) is universal at s — > oo. In particular, 

/||(x) = /^a -1 / 2 (s) yjj(xi)^ • ^a _1 / 2 (s) fi(xi + x)^ \ with i ^ 1 is a universal function. 
According to we have Gy(x) = 5(s)x A ll~ d = (n - l)(2M)" 2 a 2 (s) [/||(x) - /||(oo) 



at 



x 



oo and s — > oo. The universality of /y(x) then implies that b{s)M ja (s) must be uni- 
versal at s — > oo. According to the scaling rules a(s) = a(l) • s x± ~ d and 5(s) = 5(1) • s A|1 
the latter is possible only if An = 2Aj_ — d holds, whence it follows also that 6(l)M 2 /a 2 (l) 
and 6(l)M 2 /a 2 (l) = bM 2 /a 2 are universal, i.e., dependent merely on n and d. Thus we 
recover one of relations (|62|). as well as the universality of the ratio bM 2 /a 2 discussed in 
Sec. IP 

It is very likely that the accuracy of (|78|) is limited even at s — > oo. However, there 
exists a less constrained form 

T s (H/T) =Y,Q{(M- </?i(x)) Ma-^s), a~ 1 /2 (s ) ^( x )| (81) 

X 

for the renormalized Hamiltonian at s — > oo, as consistent with the idea that Eq. (|79|) 
with (po = M + ((fj_/M) holds approximately for relevant configurations of block-spins. 
It also leads to the relation An = 2Aj_ — d and the universality of bM 2 /a 2 . 

Contrary to the previous discussion in Sec. 19.21 in this case our conclusions are not 
restricted to small u. 
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9.6 Magnetization in a small external field 

Here we discuss the qualitative behavior of magnetization M in a small external field. 

The modulus conservation principle holds far below the critical point (at large negative 
ro ) where the fluctuations of modulus | y?(x) | are reduced to a small vicinity of ipo(h) — 
\J— ro/(2n) — H/{4tq), as consistent with the minimum of Hamiltonian In this case 
we have 

M{h) = <^(x)) ~ Mh) (^ 2 (x)) (fc)) , (82) 

where 0(x) is the angular deviation of y(x) from the direction of the external field h. We 
consider the limit lim lim lim where Eq. (1821) is asymptotically exact, since I </?(x) l = 

c—>oo r ^-oo h— »0 

const holds with an unlimited accuracy and, simultaneously, the angular fluctuations are 
suppressed. The variation of (fo(h) with h is analytical, whereas the singular behavior of 
M{h) at h — > is due to the term 

(fl 2 (x)) (h) ~ Gj_(0) = (2vr)- d J G ± (k) dk . (83) 

The transverse correlation function behaves like Gj_(k) ~ afc~ A± (the k — ► asymptotic 
at h = 0) when k is decreased below some k cr u (defined for any given ro and c) if h — > 
until it saturates at the known value M/h valid for k = 0. From this we find 

M(h) - Af (+0) oc fc^A-L)- 1 at /i -> . (84) 

Since the exponent p = (d/X±) — 1 is universal, Eq. (|84|) is valid for any T < T c including 
vicinity of the critical point. This yields the longitudinal susceptibility 

X\\ = dM(h)/dh oc /i( d -2*x)M± = h~ x \\/ x± at h->0. (85) 

According to (|62j) . we have (d/2)— 1 < p < 1, which yields 1/2 < p < 1 in three dimensions. 
The lower value p = 0.5 corresponds to the conventional statement |13 | ll7 | l24"| that Aj_ = 2. 
Our numerical test in the following section, however, does not support this possibility. 



10 Monte Carlo test in 3D XY model 

To verify the theoretical predictions for the exponent p = (d/X±) — 1 in l|84jl. we have 
made Monte Carlo (MC) simulations of 3D XY model on simple cubic lattice with the 
Hamiltonian 

l = -^(E s ^ + E hs i] > ( 86 ) 

where Sj is the spin variable (two-component vector) of the i-ih lattice site, K is the cou- 
pling constant, and h is the external field. Based on the universality argument (cf. Sec. l9.5|) . 
the same value of p is valid also for the actual <p model with two-dimensional order 
parameter (n = 2). The simulations have been made in the ordered phase at K = 
0.475,0.5,0.55 > K c , where K c ~ 0.4542 is the critical point. Only the case K = 0.5 
is discussed in detail, since the analysis made at K = 0.475 and K = 0.55 is similar. 
The quatity (| m |), where m is the magnetization per spin, has been evaluated for dif- 
ferent linear sizes of the lattice L. The Wolff's cluster algorithm has been used at 
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Table 1: The MC simulated values of (| m |) for 3D XY model depending on the external field at 
a fixed coupling constant K = 0.5. 



11 


(\m\) 


L=8 


L=16 


L=32 


L=64 


0.028 


0.586155(85) 


0.567394(63) 


0.562201(78) 


0.561218(45) 


0.04 


0.593381(86) 


0.576368(91) 


0.572118(70) 


0.571425(55) 


0.056 


0.601719(100) 


0.586635(99) 


0.583529(77) 


0.583020(45) 


0.08 


0.613203(93) 


0.600377(109) 


0.597872(63) 


0.597552(52) 


0.112 


0.626239(89) 


0.615581(81) 


0.613741(61) 


0.613527(46) 


0.16 


0.643071(78) 


0.634658(63) 


0.633304(50) 


0.633242(31) 


0.224 


0.661554(58) 


0.655097(58) 


0.654166(34) 


0.654016(28) 


0.32 


0.684105(62) 


0.679169(51) 


0.678533(31) 


0.678483(22) 


0.448 


0.707654(44) 


0.704018(37) 


0.703606(32) 


0.703496(20) 


0.64 


0.734667(38) 


0.732058(25) 


0.731754(22) 


0.731714(16) 



h = 0. The results are (| m |) = 0.570297(23), 0.542411(21), 0.530317(20), 0.524606(19), 
0.521846(26), and 0.520449(35) for L = 8, 16, 32, 64, 128, and 256, respectively. These 
values for L = 8 to 256 have been obtained by an averaging over, respectively, 8 • 10 7 , 
2 • 10 7 , 5 • 10 6 , 1.25 • 10 6 , 3.125 • 10 5 , and 7.5 • 10 4 cluster algorithm steps. Each simulation 
has been split typically in 51 bins to calculate the mean value and the standard deviation 
discarding the first bin (first 2 bins at L = 256). The value at L = 256 has been obtained 
by a weighted averaging over two simulations including totally 30 (10 + 20) not discarded 
bins, each consisting of 2500 cluster updates. 

The simulations at h > have been done by the Metropolis algorithm. The results 
of simulation for L = 8, 16, 32 and 64 are listed in Tab. ^ The statistical averages have 
been evaluated from 3.2 • 10 6 x (64/L 2 ) sweeps at 0.08 < h < 0.64, discarding no less than 
50 000 sweeps from the beginning of the simulation to ensure an accurate equilibration. 
The total length of the simulation as well as the discarded part have been increased by 
a factor 0.08/Zi at h < 0.08. Like in the case of h = 0, each simulation has been split in 
bins, using the last 50 ones for the estimations. 

The linear congruatial generator with multiplier 7 5 and modulo 2 31 — 1 (Lewis genera- 
tor), improved by a shuffling, has been used as a source of pseudo random numbers. The 
standard shuffling scheme jHS] with the length of string N = 10 6 has been completed by 
a second shuffling, where the whole cycle (consisting of 2 31 — 2 numbers) of the original 
generator has been split in 2 20 segments, restarting the generation from the beginning of a 
new randomly choosen segment when the previous one is exhausted. The scheme provided 
excellent results in test simulations of 2D Ising model, where the simulated values of in- 
ternal energy, specific heat Cy, and its first two derivatives have been compared with the 
exact results. No systematic deviations have been observed in rather long simulations at 
the critical point providing the standard error in Cy as small as 0.02% for 48 x 48 lattice 
and 0.11% for 256 x 256 lattice. 

The quantity M(+0) in has been evaluated by extrapolating our (| m |) data to 
the thermodynamic limit, based on an empirical observation that (| m |) is almost linear 
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Figure 1: The mean magnetization modulus (\m\) vs the inverse linear system size 1/L at 
K = 0.5, h = 0. The statistical errors are much smaller than the symbol size. The linear 
extrapolation (dashed line) and the quadratic fit (solid line) yield the asymptotic values 
0.519073(37) and 0.519116(30), respectively. 



function of the inverse lattice size 1/L, as it is evident from Fig. ^ The linear fit of 
three largest sizes (dashed line in Fig. provides M(+0) = 0.519073(37), the quadratic 
fit within 16 < L < 256 yields M(+0) = 0.519116(30) (solid line), whereas the cubic fit 
within 8 < L < 256 gives us M(+0) = 0.519096(33). The values of goodness Q (see [HZI 
for the definition) of the respective fits are 0.769, 0.767, and 0.794. We have accepted the 
result M(+0) = 0.519116(30) of the quadratic fit for our further estimations, although 
we could choose any of two other values as well, since the differences are rather small. 
Moreover, even a shift by 20 standard deviations would not change the qualitative picture 
in Fig. [2 where the effective critical exponent p e ff(h,L) is shown, defined as the mean 
slope of the hx[M(h',L) - M(+0)] vs ln/i' plot within h! G [h,2h], where M(h',L) is the 
value of (| m |) at the field b! and the linear lattice size L. Thus, the effective exponent is 
given by 

~M(2h,L) -M(+0)' 



p e ff(h,L) = In 



/In 2, (87) 



M(h,L) - M(+0) 

and the true value of the critical exponent p is obtained in the limit p = lim lim p e tt{h, L). 

h— »0 L— >oo 

As we see from Fig. [21 the effective exponent p e ff(h, L) increases monotoneously and con- 
verges to a certain value p e ff(h) = lim p e ff(h, L) with increasing of L. The data in Tab.^ 

obey an approximate scaling relation p e ff(h,2L) — p e ff(h,L) = 
A[p e ff(4h,L) — p e ff(4h, L/2)] which holds with an almost constant value of coefficient 
A, i. e., A ~ 1.15 for L = 16 and A ~ 1.2 for L = 32. According to the finite-size 
scaling theory L/£ is the only essential scaling argument. It implies that the correlation 
length ^ is roug hly proportional to /i~ 1/2 within h G [0.028; 0.64]. At h = 0.112 our results 
for L = 32 closely agree with those at L = 64, which indicates that £ is several times 
smaller than 32 in this case. The actual scaling analysis then implies that £ is several 
times smaller than L = 64 at h = 0.028. The scaling relation we found allows us to 
evaluate p e ff(h, 64) — p e fj(h) for h > 0.028. It leads to the conclusion that, with a high 
enough accuracy, our results at L = 64 (solid circles in Fig. EJ) correspond already to the 
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Figure 2: The effective exponent p e ff(h,L) for K = 0.5 evaluated at L = 16 (triangles), 
L = 32 (empty circles) , and L = 64 (solid circles) . The squares and rhombs show the L = 
64 results for two other couplings (temperatures) K = 0.475 and K = 0.55, respectively. 
The estimates for L = oo, obtained by correcting the L = 64 values, are swown by pluses. 
The statistical errors are within the symbol size. The dashed and solid lines show the fits 
p e ff(h) = p + aih 1 / 2 + d2h and p e ff(h) = p+ah u , respectively. The horizontal dot-dashed 
line shows the asymptotic value of the "exact" RG theory. 

thermodynamic limit at not too small fields h > 0.04, whereas a small correction (incre- 
ment) about 0.002 is necessary to get the true value of p e ff(h) at h = 0.028. This value is 
indicated in Fig. [21 by a plus. Even larger than L = 64 lattices have to be considered for 
a reliable estimation of the thermodynamic limit at h = 0.02 an smaller fields. 

According to the "exact" RG theory critised in Sec. 19.11 one can expect that p e tt{h) 
converges to the asymptotic value p = 0.5 linearly in h 1 / 2 at h — ► 0, as consistent with 
the expansion of M in powers of h 1 ^ 2 . We have used this scale in Fig. [21 to show that 
the p e ff(h) vs h 1 / 2 plot at K = 0.5 (solid circles at h > 0.28 and plus at h = 0.28) goes, 
indeed, almost linearly, but clearly not to the "exact" value 0.5. Two different fits we 
made, namely, p e ff(h) = p + a\h 1 / 2 + 02/1 (dashed curve) and p e ff(h) = p + ah w (solid 
curve) yield p = 0.705(9) and p = 0.728(25) (with uj = 0.394(50)), respectively. The fits of 
the first kind have been made also at K = 0.475 and K = 0.55 (lower and upper dashed 
lines) providing p = 0.725(9) and p = 0.707(36), respectively. The fits of the second kind 
are not stable enough in these cases because of too large inaccuracy in u. The above listed 
estimates of p satisfactory well agree with the average value about 0.716 and, thus, confirm 
the expected universality of this exponent. A similar method of effective exponents has 
been tested in 2D Ising model ^2] where it provided very accurate results. 

The actual results for XY model agree with our prediction 1/2 < p < 1 for 3D case 
(cf. Sec. 19. 6|) and are incompatible with the conventional believe that p should be 1/2. 
To the contrary, it has been claimed in |38| I39j that the MC simulated magnetization 
data for 0(4) and 0(2) models well agree with the predictions of the Gaussian spin wave 
theory. However, we failed to see any serious argument in these papers by J. Engels et. 
al, since the only quantity which could be extracted from the magnetization data and 
precisely compared to the theory, i. e., the universal exponent p, has not been evaluated 
there. Moreover, their magnetization plots for the 0(2) model are remarkably nonlinear 
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functions of h 1 / 2 , i.e., they do not provide even an indirect evidence that p is just 1/2. The 
only values we found in [HHj, which can be compared to ours, are M(+0) = 0.5186(01) 
at K = 0.5 and M(+0) = 0.6303(1) at K = 0.55. Our respective values 0.519116(30) 
and 0.630545(24) are similar. Since there is no reason to assume that one of us has made 
a wrong simulation, the small but remarkable (5 their standard deviations at K = 0.5) 
discrepancies, obviously, are due to the different extrapolation procedures used. Our values 
are more precise and reliable, since our method allows to estimate M(+0) from simulations 
at h = without any extrapolation over h used in [35]. Besides, we have used a larger 
maximal system size L = 256 as compared to L = 160 in [22] . As a result, the extrapolation 
gap in M is by an order of magnitude smaller in our case. The spontaneous magnetization 
M(+0) has been evaluated in jHHl from the fit M(h) = M(+0) + ah p + bh with p 

I 2. In such away, the extrapolation which is biased by p = 1/2 gives systematically 
underestimated values of M(+0) thus providing an indirect evidence that p / 1/2. If 
p = 1/2 is replaced with p > 1/2 in this ansatz, then the extrapolation gap becomes 
smaller (M(+0) becomes larger) and the above discussed small discrepancies in M (+0) 
can be removed. This analysis suggests that p > 1/2 holds, as consistent with the direct 
estimation in our paper. 

II Conclusions 

In the present work we have extended our diagrammatic method introduced in .5] to study 
the ip A model in the ordered phase below the critical point, i. e., at T < T c . In summary, 
we conclude the following. 

1. The diagrammatic equations derived in [£] have been generalized to include the 
symmetry breaking term fixing the axis of ordering at T < T c (Sec. |2J). 

2. An alternative formulation of our equations has been proposed. It has been shown 
that our equations coincide with the free energy variation principle ( Sec. 01) ■ 

3. The solution for the two-point correlation (Greens) function depending on tem- 
perature T has been analyzed qualitatively not cutting the perturbation series. It 
includes the low-temperature solution at tq — ► — oo in the case of scalar order- 
parameter field (Sec. IIJ), as well as the general solution of n-component vector model 
at T -> T c (Sec. El). 

4. Based on our diagrammatic equations, the asymptotic long-wave (k — » 0) behavior 
of the transverse and longitudinal Greens functions below T c has been analyzed in 
Sees. El to El This analysis shows that G±(k) ~ ak~ x± and Gy(k) ~ bk~ x W 
with exponents d/2 < Aj_ < 2 and Ay = 2Aj_ — d, and with universal ratio bM 2 Jo? 
is the physical solution of our equations at the spatial dimensionality 2 < d < 4, 
which coincides with the asymptotic solution at T — > T c as well as with the known 
rigorous results discussed in Sec. EUand the non-pertubative renormalization group 
arguments provided in Sec. 19.51 It is confirmed also by the Monte Carlo simulations 
in Sec. 1101 Formally, the results of the perturbative RG theory are recovered at 
Aj_ = 2. However, we have disproven the conventional statement of this theory that 
Aj_ = 2 is the exact result (Sec. 19.11) . 
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5. Monte Carlo simulations in 3D XY model has been performed (Sec. I1U|) to test the 
exponent p, describing the behavior of the magnetization M(h) — M(+0) ~ h p in 
small external fields h — ► below T c . The simulation results confirm the universality 
of this exponent as well as our prediction 1/2 < p < 1, and are incompatible with 
the value p = 1/2 of the Gaussian spin wave theory supported by the perturbative 
RG analysis. 



Appendix 

Here we study the thermodynamic limit of E(k, £), starting with k = 0. For simplicity, 
we consider only the case of scalar order-parameter-field, i. e., n = 1. The extension to 
the n-component case is trivial: the fixed zero-vectors always refer to the longitudinal 
component. According to the definition of S(0, C) (cf. Eqs. (jllj) to (|T1)0 . this quantity 
obeys a selfconsistent diagrammatic equation 



e(o,o :! -°-0- - + -KS>- - 



(-41) 



where the block 



represents a resummed perturbation series of all connected 



diagrams of this kind, which do not contain parts like ^ and/or v/xo^xru and 

1 



cannot be split in two as follows 



. At C 
k 4 = -ki 



with four fixed outer lines and vectors ki, k2, k3 
perturbation sum of the cumulant four-point correlation function 

G c (k 1 ,k 2 ,k 3 ) = G(k 1 ,k 2 ,k 3 ) 



-G (ki, k 2 ) G (k 3 , k 4 ) - G (k a , k 3 ) G (k 2 , k 4 ) - G (k 1} k 4 ) G (k 2 , k 3 



all connected diagrams 
k 2 — k 3 represent the 



(A2) 



where G(ki,k 2 ) = <5k x _k 2 G ! (ki) is the two-point correlation function. In principle, the 
diagrams of G c (ki, k 2 , k 3 ) can be grouped like those of S(k, Q following It implies the 
summation over the chains of blocks — ^ — , yielding skeleton diagrams with respect 

to the solid lines, followed by the summation over the chains of blocks , yielding 

skeleton diagrams with respect to the waved lines. We remind that the skeleton diagrams 
are those connected diagrams where the true correlation function G(k) is related to the 
solid lines and the dashed lines inside the diagrams are replaced by the waved lines. 
Besides, the skeleton diagrams do not contain parts — ^ — and >a/^w . In such 

wyy 

a way, the perturbation sum of G c (0,0,0) is almost the same as n KxSj n with the 



only difference that it includes also the diagrams like 
the equation for G c (0, 0, 0), 



G c (0,0,0) 



o o- 

.0 + 0- 



X 



It leads to 



(A3) 



in which the diagram blocks are calculated at £ = 1. In analogy to G(0) ~ M 2 V 
(cf. eq. (|23|0 . the zero- vector term V~ 2 G (0,0,0) represents the constant (non-decaying) 
contribution M to the real-space four-point correlation function below T c . According 
to (A2), it means that 



G c (0, 0, 0) = -2M 4 V 2 at V 



oo 



(A4) 
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At the first step we take into account only those contributions where k = vectors are 
related to the connecting solid lines inside the last diagram in (A3). Then, together 
with (A4), (Al), and (|12|). we have a system of selfconsistent equations yielding 

£(0,1) ocV (A5) 

with two possible values of the proportionality factor — ^5 ± Vvfj uM 4 . One of possi- 
bilities is that the additional terms with non-zero wave vectors related to the connecting 
solid lines in (A3) give a negligible corretion at V — > oo, as consistent with the following 
idea. The four-point function G(k, — k, 0) behaves like M 2 FG(k) below T c , as it follows 
from a simple consideration of the limit where two points of the real-space four-point 
correlation function are infinitely distant. However, this term cancellates in the cumulant 
average G c (k, — k, 0) at k ^ entering the equation, analogous to (A3), for the block 



-k i88 o • As a result, a selfconsistent solution exists where this block is vanishingly 

small as compared to V and Eq. (A5) remains correct. Nevertheless, other kind of selfcon- 
sistent solutions cannot be excluded where all terms on the right-hand side of Eq. (A3) 



k- 



and of similar equations for the blocks -k kxSI n and -k kx>l q are compatible 



■0 k- 



■q 



The latter is possible if the block in (A3), having four fixed zero-vectors, diverges as V 2+ ^ 
with /x > 0, whereas the above mentioned blocks with non-zero vectors k and q diverge 
like V 1+>1 and V^, respectively. In this case S(0, 1) is proportional to V 1+ ^ . 

Contrary to the case of k = 0, the thermodynamic limit exists for the function S(k, 1) 
at k ^ 0, as it can be found easily by a suitable grouping of the divergent terms in 
Only those terms can diverge at V — > oo which contain insertions with 2m outer solid 



y/x o_K7/j_ () 

zero- vector lines, like n v/X n , o^^^o j etc. In an ordinary case m factors 



G(0) ~ M 2 V related to m solid lines with fixed k = vectors are compensated by 
a removal of m integrations over wave vectors, as consistent with the well known rule 
J2k —* V(2ir)~ d J dk. This condition is violated if the constraints k = are not indepen- 
dent. The sum of all wave vectors coming out from any node is zero. As a consequence, 
this property holds also for any block. Therefore only 2m — 1 constraints k = for 2m 
outer solid lines of the above discussed insertions are independent, i. e., only 2m — 1 in- 
tegrations are removed. As a result, any such insertion provides a diverging factor oc V 
for the resulting diagram unless this factor is compensated by vanishing waved line with 
fixed k = 0. Due to these properties, the constraints k = for 2m solid lines connecting 

two parts of a diagram like ^ also are not independent, but this situation 
is possible only for the diagrams of S(0, £)• 

Thus, some diagrams containing specific insertions with 2m outer zero-vector solid 
lines are divergent. At the same time, all such insertions with 2m outer lines represent the 
perturbation sum of the 2m-point cumulant correlation function with all zero arguments. 
In general, the 2m-point cumulant G^ m ^ (ki, k 2 , • • • , \^2m-i) is defined as 

Gi 2m ) (ki , k 2 , . . . , k 2m _x ) = (ki , k 2 , . . . , k 2m _i ) - (ki , k 2 , . . . , k 2m _! ) , (A6) 

where G^ m ^ is given by resummed connected diagrams, i. e., a diagram block with 2m 
fixed outer solid lines, whereas S^ 2 " 1 ) represents the sum over all possible splitings of this 
block in smaller parts. These functions contain only 2m — 1 independent arguments (wave 
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vectors of outer solid lines), since the sum over all 2m wave vectors is zero. The correlation 
function G^ m "> (0, 0, ... , 0) is related to the non-decaying part M 2m of the real-space 2m- 
point correlation function and, thus (according to the Fourier transformation), is M 2m V m 

at V —> oo. From (A6) we find that G^™^ (0,0,..., 0) is proportional to M 2m V m . In such 
a way, we can resummate the specific zero-vector insertions to replace the perturbation 
sums by corresponding 2m-point cumulants. It is then straightforward to see that any 
inserted 2m-point zero- vector cumulant does not cause a divergence at V — ► oo, but only 
renormalize the original diagram by a constant factor. Namely, m solid zero- vector lines of 
the original diagram are broken to insert the cumulant block, thus replacing the previous 
factor [G(0)] m ~ M 2m V m of these zero-vector lines with the cumulant value Q m M 2m V m , 
where Q m is a constant. 

Not only the function S(k, Q, but also the derivative <9X(q, ()/dG(k) is relevant to our 
analysis. This derivative is represented by diagrams obtained by breaking one solid line 
with wave vector k in the diagrams of S(q, £), removing the corresponding factor G(k). 
Let us consider those diagrams of 9S(q, l)/dG(k) with k ^ 0, containing the zero-vector- 
cumulant insertions, where the broken solid line does not belong to any of the inserted 
blocks. The diagrams including resummed blocks of this kind and the corresponding 
"normal" diagrams whith no insertions differ merely by constant factors. Namely, as in 
the case of S(k, 1), any insertion of a resummed 2m-point cumulant block with 2m outer 
zero-vector lines gives a constant factor Q m at V — ► oo. 

Below we prove that the derivative dQ m /dG(k) vanishes at V — > oo. As a result, the 
corresponding terms where solid line is broken inside a zero-vector-cumulant block do not 
give an extra contribution to 9S(q, l)/9G(k). Quantity M 2m V m dQ m / dG(\t), is repre- 
sented by the skeleton diagrams of G^ m ' 1 (0,0,..., 0) in which one solid line with vector k 
is broken in two, removing the factor G(k). If factors G(k) and G(— k) = G(k) are restored 
for both parts of the broken line, we obtain the diagrams 

of G (2m+2) ( kj _ kj0> 0) To 

simplify the further notation we shall replace the above set of arguments with one argu- 
ment k. Thus, a resummation of these diagrams yields 

M 2m V m [G(k)] 2 dQ m /dG{k) = Gf m+2) (k) . {AT) 

According to the physical arguments we have 

G (2m) {0) = M 2m V m , (AS) 

G (2m+2) = g qq M 2mym ^g) 

at V — > oo. Eq. (A8) represents the condition for the non-decaying part M 2m V m of 
the 2m-point real-space correlation function, whereas (A9) describes the limit where 2m 
points of the (2m + 2)-point function are infinitely distant. By setting m = 2 in (A6) and 
taking into account (A9), we obtain at k / 

M-V^G^k) = G(k) - G(k) G(0)M- 2 r' = at V ->■ oo . (AW) 

We can prove by induction over £ that M~ 2e V~ e G^ e+2 \k) = holds at V —> oo for any 
£ > 1 and k / 0. Eq. (A10) means that it holds at £ = 1. If it holds for £ < m, then the 
only relevant terms in the equation (A6) for the (2m + 2)-point cumulant are 

G i2m+2) (k) = G (2m+2) (k) - G(k) G i2m) (0) - G(k) S {2m) (0) , (All) 
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which yield vanishing result for M~ 2m V~ m G c m (k) according to Eq. (A6) for zero- 
vector cumulants, as well as (A8) and (A9). According to (A7), the latter means that 
quantity dQ m /dG(k.) vanishes (at k 7^ 0) in the thermodynamic limit V — > 00. 

However, if it would not vanish, then it would produce a contribution to the derivative 
<9£(q, l)/dG(k) which is by a factor V larger than the ordinary terms ~ V^ 1 . It means 
that not only the leading behavior, but also the corrections to finite-size scaling in (A8) 
and (A9) could play some role. At n = 1 these corrections are exponentially small, as 
consistent with the known exponential decay of the real-space correlation functions in 
Ising model below T c . The latter is true also at n > 1 for the longitudinal component of 
the correlation functions at finite values of the amplitude 5 of the symmetry-breaking term 
in (j2J), since in this case Hamiltonian (j2J) looses the rotational symmetry and, therefore, 
belongs to the Ising universality class. Since 5 is a continuous parameter, the exponential 
decay of correlations will be observed at large enough distances x for any arbitrarily small, 
but finite value of 5. It means that the corrections are exponentially small in our limit 

lim lim , where L is the linear size of the system. There is no any contradiction with the 

S— >0 L^oo 

power-like decay of correlations we found at n > 1, since this decay refers to the limit 
where x — ► 00 and, simultaneously, x/L — > hold (i. e., we find the thermodynamic limit 
at any given nonzero wave vector). 

Thus, we finally arrive to the conclusion that corrections to (A8) and (A9) are ex- 
ponentially small in L = V l / d and, therefore, dQ m /dG(k) = holds at k / with a 
high enough accuracy, i. e., quantities Q m can be treated as pure constants providing no 
extra contributions due to their variations. As discussed before, it means that the in- 
serted zero-vector-cumulants merely renormalize by constant factors the "normal" terms 
in the diagram expansion of S(k, 1) and 5S(q, l)/<9G(k). Since our analysis is not based 
on specific values of expansion coefficients, this renormalization cannot affect our quali- 
tative conclusions. It means that the general scaling form of E(k, 1) and <9£(q, l)/<9G(k) 
is correctly predicted by the simplified analysis which ignores the zero-vector-cumulant 
insertions. 

Our foregoing consideration of S(k, £) is restricted to £ = 1 due to a technical reason 
that the relation to cumulant correlation functions with a certain physical meaning is 
known only for this case. Since £ is a continuous parameter, we believe that the above 
discussed properties remain true for all values of interest, i. e., ( G [0; 1]. 

We can use the alternative method to find quantity i?j(k), as proposed in Sec. 01 It 
is straightforward to check in each specific case we considered in our paper that Eq. ifTTjl 
provides the same scaling properties of Ri(k) at k 7^ as the equations (jHJ) and (j3Ul> . 
Regarding the "normal" terms, it is easy to verify (like in Sees. 19.21 - 19.31 and in Sec. EJ) 
that in the most nontrivial case of n > 1 at k — > below T c , as well as at ro — * —00 in 
the case of n = 1, the main contribution comes from all diagrams of Eq. (|17|) in which 
the waved lines have nonzero wave vectors and each second solid line with i = 1 has fixed 
wave vector k = (i. e., there are no integrations in the case of n = 1), counting the 
pair of outer lines as one line. This peculiarity does not refer to the case T — * T c , where 
all the partial contributions with M°, M 2 , M 4 , etc., are compatible. The two different 
modifications of our method obviously give consistent corrections to scaling: they have the 
same origin. The alternative approach [Eq. (|17|) ] does not suffer from the problems with 
parameter £ < 1, since we have £ = 1. As before, the zero-vector-cumulant insertions 
merely renormalize by constant factors the "normal" terms. It proves the statement that 
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the simplified analysis, ignoring these insertions, provides correct general scaling form of 
the solution for Ri(k) and, therefore, Gj(k). 
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